home *** CD-ROM | disk | FTP | other *** search
/ FishMarket 1.0 / FishMarket v1.0.iso / fishies / 051-075 / disk_071 / airfoil / airfoil.c < prev    next >
C/C++ Source or Header  |  1992-05-06  |  13KB  |  484 lines
  1. /********************************************************************
  2. *                                                                   *
  3. * Joukowski Airfoil Generator with Streamline and Pressure          *
  4. * Distribution Algorithms                                           *
  5. *                                                                   *
  6. * Written by:   Russell Leighton      15 March 1987                 *
  7. *               Lancaster, CA                                       *
  8. *                                                                   *
  9. ********************************************************************/
  10.  
  11. #include "airfoil.h"
  12.  
  13. main()
  14. {
  15.    float rs,theta,h,vel;
  16.    int psi;
  17.    ULONG MessageClass;
  18.  
  19.    open_things();
  20.    do_about();
  21.  
  22.    /****************/
  23.    /* Loop forever */
  24.    /****************/
  25.  
  26.    for(;;)
  27.    {
  28.       while (Continue)
  29.       {
  30.          /****************************************************/
  31.          /* Wait, initially and after each plot, for user to */
  32.          /* bring up the double menu requester               */
  33.          /****************************************************/
  34.  
  35.          Wait(1L<<w->UserPort->mp_SigBit);
  36.          if((message = (struct IntuiMessage *)GetMsg(w->UserPort)) != NULL)
  37.          {
  38.             MessageClass = message->Class;
  39.             ReplyMsg(message);
  40.             if (MessageClass == REQVERIFY)
  41.             {
  42.                do_request();
  43.                break;
  44.             }
  45.          } /* end if */
  46.       } /* end while */
  47.  
  48.       do_init();
  49.  
  50.       if(mode)
  51.       {
  52.  
  53.          /********************/
  54.          /* Plot Streamlines */
  55.          /********************/
  56.  
  57.          FILL = TRUE;
  58.          for (psi = 12;psi > 0;--psi)
  59.          {
  60.             do_mess();
  61.             if(!Continue) break;
  62.  
  63.             PENUP;
  64.             SetAPen(rp,(long)(psi+1));
  65.             SetBPen(rp,(long)(psi+1));
  66.  
  67.             for (theta = 0.015;theta < TWO_PI;theta += PI/100)
  68.             {
  69.                vel = 2*velocity*sin(theta);
  70.                vel = abs(vel);
  71.                rs = (psi+sqrt(psi*psi+r*r*vel*vel))/vel;
  72.                transform(rs,theta);
  73.             } /* end for */
  74.             AreaEnd(rp);
  75.          } /* end for */
  76.  
  77.          /*******************************/
  78.          /* Plot Stagnation Streamlines */
  79.          /*******************************/
  80.  
  81.          do_mess();
  82.          if(Continue)
  83.          {
  84.             FILL = FALSE;
  85.             PENUP;
  86.             SetAPen(rp,1L);
  87.             h = (r-4.0)/40.0;
  88.             theta = 0.0;
  89.  
  90.             for (rs = 4;rs >= r;rs += h)
  91.             {
  92.                transform(rs,theta);
  93.             }
  94.  
  95.             PENUP;
  96.             theta = PI;
  97.  
  98.             for (rs = 4;rs > r;rs += h)
  99.             {
  100.                transform(rs,theta);
  101.             }
  102.          } /* end if */
  103.       } /* end if */
  104.  
  105.       else
  106.       {
  107.  
  108.          /******************************/
  109.          /* Plot Pressure Distribution */
  110.          /******************************/
  111.  
  112.          do_mess();
  113.          if(Continue)
  114.          {
  115.             FILL = TRUE;
  116.             PENUP;
  117.             SetAPen(rp,2L);
  118.             SetBPen(rp,2L);
  119.  
  120.             for (theta = 0.0;theta <= TWO_PI;theta += PI/100)
  121.             {
  122.                rs = r+sin(theta)*sin(theta);
  123.                transform(rs,theta);
  124.             } /* end for */
  125.             AreaEnd(rp);
  126.          } /* end if */
  127.       } /* end else */
  128.  
  129.       /****************/
  130.       /* Plot Airfoil */
  131.       /****************/
  132.  
  133.       do_mess();
  134.       if(Continue)
  135.       {
  136.          FILL = TRUE;
  137.          PENUP;
  138.          rs = r;
  139.          SetAPen(rp,1L);
  140.          SetBPen(rp,1L);
  141.  
  142.          for (theta = 0.0;theta <= TWO_PI;theta += PI/100)
  143.          {
  144.             transform(rs,theta);
  145.          }
  146.          AreaEnd(rp);
  147.       } /* end if */
  148.    } /* end forever */
  149. } /* end main */
  150.  
  151. do_init()
  152. {
  153.    float a0;
  154.  
  155.    SetAPen(rp,0L);
  156.    RectFill(rp,(long)(XMIN+1),(long)(YMIN+1),(long)(XMAX-1),(long)(YMAX-1));
  157.    SetOPen(rp,1L);
  158.  
  159.    /***********************************************************/
  160.    /* Calculate circle constants (circle center and radius)   */
  161.    /* from airfoil constants through a reverse transformation */
  162.    /***********************************************************/
  163.  
  164.    alpha = (float)angle*PI/180;
  165.    c = (float)camber/25;
  166.    t = (float)thickness/25;
  167.    a = 0;
  168.    b = c/2;
  169.  
  170.    do
  171.    {
  172.       a0 = a;
  173.       a = t*(2*a0+1)/4/sqrt(b*b+2*a0+1);
  174.       b = c*(1+2*a)/(2+2*a);
  175.    }
  176.    while(abs(a-a0) > TOL);
  177.    
  178.    r = sqrt(b*b+(a+1)*(a+1));
  179.    Continue = TRUE;
  180. }
  181.  
  182. do_mess()
  183. {
  184.    ULONG MessageClass;
  185.  
  186.    /******************************************************************/
  187.    /* Check for double menu requester. This can be brought up at any */
  188.    /* time but can not be displayed unless the message is replied    */
  189.    /* to, therefore this routine must be called periodically.        */
  190.    /******************************************************************/
  191.  
  192.    if((message = (struct IntuiMessage *)GetMsg(w->UserPort)) != NULL)
  193.    {
  194.       MessageClass = message->Class;
  195.       ReplyMsg(message);
  196.       if(MessageClass == REQVERIFY) do_request();
  197.    }
  198.    ixo = XCEN;
  199.    iyo = YCEN;
  200. }
  201.  
  202. transform(rs,theta)
  203.  
  204. float rs,theta;
  205. {
  206.    float x,y,z,u,v;
  207.    int PLOT;
  208.    long ix,iy,cx,cy;
  209.  
  210.    /********************************************************************/
  211.    /* This is the Joukowski transformation routine. This is also where */
  212.    /* the plotting is done. Plotting is usually done using the area    */
  213.    /* fill routines (AreaDraw & AreaMove). If FILL is true then the    */
  214.    /* points are used to build up the area shape to be filled. See the */
  215.    /* Rom Kernal manual containing the graphics primatives for more    */
  216.    /* details. If FILL is false then normal plotting is done.          */
  217.    /********************************************************************/
  218.  
  219.    x = rs*cos(theta-alpha)+a;
  220.    y = rs*sin(theta-alpha)+b;
  221.    z = 1/(x*x+y*y);
  222.    u = x*(1+z)*cos(alpha)-y*(1-z)*sin(alpha);
  223.    v = y*(1-z)*cos(alpha)+x*(1+z)*sin(alpha);
  224.  
  225.    ix = (long)(SFAC*u+XCEN);
  226.    iy = (long)(YCEN-SFAC*v);
  227.  
  228.    if(FILL)
  229.    {
  230.       PLOT = FALSE;
  231.       cx = ix;
  232.       cy = iy;
  233.  
  234.       /*************************************************************/
  235.       /* This subroutine also clips the display area, hence the    */
  236.       /* following statements. Contrary to the what the Rom Kernal */
  237.       /* manual states, all clipping must be done by the program   */
  238.       /* if the area fill routines are used otherwise a nasty      */
  239.       /* crash will result if plotting takes place out of bounds.  */
  240.       /* Only the x values are clipped accurately since, for this  */
  241.       /* routine only the x values at the boundary need to be      */
  242.       /* accurate. The PEN parameter indicates the pen status for  */
  243.       /* moves and draws as set by either PENUP or PENDOWN.        */
  244.       /*************************************************************/
  245.  
  246.       if((ix <= XMIN) && (ixo > XMIN))
  247.       {
  248.          cy += (float)(iyo-iy)*(XMIN-ix)/(ixo-ix);
  249.          cx = XMIN;
  250.       }
  251.       else if((ixo <= XMIN) && (ix > XMIN))
  252.       {
  253.          iyo += (float)(iy-iyo)*(XMIN-ixo)/(ix-ixo);
  254.          ixo = XMIN;
  255.          AreaDraw(rp,ixo,iyo);
  256.          PLOT = TRUE;
  257.       }
  258.       else if((ix >= XMAX) && (ixo < XMAX))
  259.       {
  260.          cy += (float)(iyo-iy)*(XMAX-ix)/(ixo-ix);
  261.          cx = XMAX;
  262.       }
  263.       else if((ixo >= XMAX) && (ix < XMAX))
  264.       {
  265.          iyo += (float)(iy-iyo)*(XMAX-ixo)/(ix-ixo);
  266.          ixo = XMAX;
  267.          AreaDraw(rp,ixo,iyo);
  268.          PLOT = TRUE;
  269.       }
  270.  
  271.       if((ixo > XMIN) && (ixo < XMAX)) PLOT = TRUE;
  272.       if(cy < YMIN) cy = YMIN;
  273.       if(cy > YMAX) cy = YMAX;
  274.  
  275.       if(PLOT)
  276.       {
  277.          if(PEN) AreaDraw(rp,cx,cy);
  278.          else { AreaMove(rp,cx,cy); PENDOWN; }
  279.       }
  280.       ixo = ix;
  281.       iyo = iy;
  282.    }
  283.    else
  284.    {
  285.       if(PEN) Draw(rp,ix,iy);
  286.       else { Move(rp,ix,iy); PENDOWN; }
  287.    }
  288. }
  289.  
  290. do_request()
  291. {
  292.    ULONG MessageClass;
  293.  
  294.    /***************************************************************/
  295.    /* This subroutine is activated when the double menu requester */
  296.    /* is brought up. It handles all input from this requester.    */
  297.    /***************************************************************/
  298.  
  299.    if(title) ShowTitle(s,FALSE);
  300.  
  301.    for (;;)
  302.    {
  303.       Wait(1L<<w->UserPort->mp_SigBit);
  304.       if((message = (struct IntuiMessage *)GetMsg(w->UserPort)) != NULL)
  305.       {
  306.          MessageClass = message->Class;
  307.          ReplyMsg(message);
  308.          switch (MessageClass)
  309.          {
  310.             case GADGETUP    :
  311.  
  312.             case GADGETDOWN  : if (do_gadgets(message) == CLOSE_GAD)
  313.                                {
  314.                                   if(title) ShowTitle(s,TRUE);
  315.                                   return;
  316.                                }
  317.                                break;
  318.          } /* end switch */
  319.       } /* end if */
  320.    } /* end forever */
  321. }
  322.  
  323. int do_gadgets(mes)
  324.  
  325. struct IntuiMessage *mes;
  326. {
  327.    struct Gadget *igad;
  328.    int gadgid;
  329.  
  330.    /*********************************************/
  331.    /* This subroutine handles all gadget input. */
  332.    /*********************************************/
  333.  
  334.    igad = (struct Gadget *)mes->IAddress;
  335.    gadgid = igad->GadgetID;
  336.  
  337.    switch(gadgid)
  338.    {
  339.       case CAMBER_GAD   : camber = (ULONG)camber_string.LongInt;
  340.                           Continue = FALSE;
  341.                           break;
  342.  
  343.       case THICK_GAD    : thickness = (ULONG)thick_string.LongInt;
  344.                           Continue = FALSE;
  345.                           break;
  346.  
  347.       case ANGLE_GAD    : angle = (ULONG)angle_string.LongInt;
  348.                           Continue = FALSE;
  349.                           break;
  350.  
  351.       case VELOCITY_GAD : velocity = (ULONG)(16-(velocity_prop.HorizPot >> 12));
  352.                           Continue = FALSE;
  353.                           break;
  354.  
  355.       case AIRFOIL_GAD  : if(mode) mode = FALSE;
  356.                           else mode = TRUE;
  357.                           Continue = FALSE;
  358.                           break;
  359.  
  360.       case TITLE_GAD    : if(title) title = FALSE;
  361.                           else title = TRUE;
  362.                           break;
  363.  
  364.       case CLOSE_GAD    : break;
  365.  
  366.       case QUIT_GAD     : close_things();
  367.                           exit(0);
  368.                           break;
  369.    }
  370.    return gadgid;
  371. }
  372.  
  373. do_about()
  374. {
  375.    int i;
  376.  
  377.    /*****************************************************/
  378.    /* This subroutine displays the initial information. */
  379.    /*****************************************************/
  380.  
  381.    SetAPen(rp,0L);
  382.    RectFill(rp,(long)(XMIN+1),(long)(YMIN+1),(long)(XMAX-1),(long)(YMAX-1));
  383.    SetAPen(rp,1L);
  384.    for(i = 0;i < 24;i++)
  385.    {
  386.       Move(rp,2L,(long)(8*i+8));
  387.       Text(rp,about[i],(long)strlen(about[i]));
  388.    }
  389. }
  390.  
  391. open_things()
  392. {
  393.    struct Library *OpenLibrary();
  394.    struct Screen *OpenScreen();
  395.    struct Window *OpenWindow();
  396.    struct ViewPort *ViewPortAddress();
  397.    BYTE *AllocRaster();
  398.  
  399.    if(!(GfxBase = (struct GfxBase *)OpenLibrary("graphics.library",0L)))
  400.    {
  401.       printf("no graphics library!!!\n");
  402.       close_things();
  403.       exit(1);
  404.    }
  405.    mask |= GRAPHICS;
  406.  
  407.    if(!(IntuitionBase = (struct IntuitionBase *)
  408.         OpenLibrary("intuition.library",0L)))
  409.    {
  410.       printf("no intuition library!!!\n");
  411.       close_things();
  412.       exit(2);
  413.    }
  414.    mask |= INTUITION;
  415.  
  416.    if (!(s = OpenScreen(&ns)))
  417.    {
  418.       printf("could not open the screen\n");
  419.       close_things();
  420.       exit(3);
  421.    }
  422.    mask |= SCREEN;
  423.  
  424.    nw.Screen = s;
  425.  
  426.    if (!(w = OpenWindow(&nw)))
  427.    {
  428.       printf("could not open the window\n");
  429.       close_things();
  430.       exit(4);
  431.    }
  432.    mask |= WINDOW;
  433.  
  434.    rp = w->RPort;
  435.    vp = ViewPortAddress(w);
  436.  
  437.    InitArea(&area,buffer,250L);
  438.    rp->AreaInfo = &area;
  439.    trp.Size = (long)RASSIZE(640L,400L);
  440.  
  441.    if (!(trp.RasPtr = (BYTE *)AllocRaster(640L,400L)))
  442.    {
  443.       printf("could not allocate memory for area fills\n");
  444.       close_things();
  445.       exit(5);
  446.    }
  447.    mask |= AREAFILL;
  448.  
  449.    rp->TmpRas = &trp;
  450.  
  451.    LoadRGB4(vp,colors,16L);
  452.  
  453.    InitRequester(&AirRequest);
  454.    AirRequest.LeftEdge = 0L;
  455.    AirRequest.TopEdge = 0L;
  456.    AirRequest.Width = 200L;
  457.    AirRequest.Height = 150L;
  458.    AirRequest.RelLeft = -100L;
  459.    AirRequest.RelTop = -75L;
  460.    AirRequest.ReqGadget = &quit_gad;
  461.    AirRequest.Flags = POINTREL;
  462.    AirRequest.BackFill = 0;
  463.    
  464.    if (!(SetDMRequest(w,&AirRequest)))
  465.    {
  466.       printf("could not set Requester\n");
  467.       close_things();
  468.       exit(6);
  469.    }
  470.    mask |= DMREQUEST;
  471.  
  472.    ShowTitle(s,FALSE);
  473. }
  474.  
  475. close_things()
  476. {
  477.    if(mask & DMREQUEST) ClearDMRequest(w,&AirRequest);
  478.    if(mask & AREAFILL)  FreeRaster(trp.RasPtr,640L,400L);
  479.    if(mask & WINDOW)    CloseWindow(w);
  480.    if(mask & SCREEN)    CloseScreen(s);
  481.    if(mask & GRAPHICS)  CloseLibrary(GfxBase);
  482.    if(mask & INTUITION) CloseLibrary(IntuitionBase);
  483. }
  484.